Optimizing the pulse shape for Schwinger pair production 
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Recent studies of the dynamically assisted Schwinger effect have shown that particle production 
is significantly enhanced by a proper choice of the electric field. We demonstrate that optimal 
control theory provides a systematic means of modifying the pulse shape in order to maximize the 
particle yield. We employ the quantum kinetic framework and derive the relevant optimal control 
equations. By means of simple examples we discuss several important issues of the optimization 
procedure such as constraints, initial conditions or scaling. By relating our findings to established 
results we demonstrate that the particle yield is systematically maximized by this procedure. 
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I. INTRODUCTION 

Creation of electron-positron pairs by a strong and 
(quasi-)static electric field, the so-called Schwinger effect, 
has been a long-standing but still unobserved prediction 
of early quantum theory [Tj-Q. All attempts to observe 
this fundamental, non-pcrturbativc effect of strong-field 
QED have failed due to its exponential suppression up 
to an unprecedented field strength of E c ~ 10 18 V/m. 
The rapid development of laser technology in recent 
years, however, has raised the hope to produce peak field 
strengths of this order such that a direct observation of 
the Schwinger effect might become possible soon 

Another possibility to create electron-positron pairs 
is to supply the necessary rest mass energy by an ar- 
bitrarily weak external field with frequency uj > 2m, 
with H = c = 1 throughout. Remarkably, this produc- 
tion mechanism works even for frequencies u> < 2m by 
an absorption of multiple photons, however, it becomes 
strongly suppressed in this regime GHz|. As a matter 
of fact, the multiphoton pair production mechanism has 
been used to explain the outcome of the SLAC E-144 
experiment more than a decade ago 

A few years ago, it was suggested to employ multi- 
photon pair production to assist the Schwinger effect, re- 
ferred to as dynamically assisted Schwinger mechanism 
Q. Note that similar proposals on how to lower the 
threshold for pair production with the aid of more ad- 
vanced setups have been put forward as well [lol - [l3j . The 
idea of the dynamically assisted Schwinger mechanism 



has been pursued recently and it has been shown that 
a proper parameter choice can in fact result in an en- 
hancement of the particle yield by orders of magnitude 
[3, [l5|. However, all these studies have been based on 
rather simple field configurations and first attempts to 
investigate more sophisticated ones have been made only 
recently [l6l - [20j . The parameter space (field strengths, 
time scales, etc.), however, grows rapidly for more com- 
plicated field configurations so that a systematic treat- 
ment becomes impracticable very soon. 

In order to overcome this drawback, we will apply op- 
timal control theory to the problem of electron-positron 
pair creation for the first time. Note that this optimiza- 
tion method has already been widely used in several other 
areas of physics such as molecular physics [2l| , quantum 
computing [22J , or splitting of Bosc-Einstcin condensates 
[23l ] , to name only a few. We will employ this method in 
order to systematically shape the electric field under cer- 
tain constraints such that the number of created electron- 
positron pairs gets maximized. In this respect, a proper 
choice of boundary conditions will be crucial. 

The structure of this paper is as follows: In section 
HIl we briefly introduce the quantum kinetic framework 
and present recent results on the dynamically assisted 
Schwinger mechanism. This section will motivate the ne- 
cessity of a systematic pulse-shaping procedure in order 
to maximize the particle yield. In section IIIII we intro- 
duce optimal control theory and apply it to the quantum 
kinetic formalism. We present the corresponding numer- 
ical results in section IIVI Finally, we conclude and give 
an outlook in section fVl 
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II. MOTIVATION 

Our electromagnetic field choice is motivated by the 
focal region of counter propagating laser pulses. Given 
that the spatial variation scale is much larger than 
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the electron's Compton wavelength all spatial inhomo- 
geneities may be ignored [24j . Moreover, assuming a 
standing wave mode we may neglect all magnetic effects 
as well, B(t) = 0. Accordingly, we will only investi- 
gate electron-positron pair creation in the presence of 
a spatially homogeneous, time-dependent electric field 
E(f) = E(t)e3. We will represent this field config- 
uration in terms of a time-dependent vector potential 
A(i) = A(t) e$ in temporal gauge A$ = 0, so that: 

E{t) = -A(i) and B(t) = V x A(t) = . (1) 
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A. Quantum kinetic formalism 

In order to calculate electron-positron pair creation in 
this type of electric background field, we will employ the 
quantum kinetic formalism. Within this approach, all 
spectral information is encoded in the distribution func- 
tion -F(q, t). We denote the canonical momentum par- 
allel to the field direction as q and restrict ourselves to 
q± = 0. Correspondingly, we will have for the transverse 
energy squared e\ = m 2 + — > m 2 . This simplifica- 
tion is justified as particle production occurs dominantly 
at small orthogonal momenta, even if a distinctive cel- 
lular structure in momentum space can be observed for 
complicated field configurations pol ]. 

As particle production is a genuine non-equilibrium 
problem in quantum field theory, it is clear that a particle 
interpretation of the distribution function is only valid at 
asymptotic times T 1 . Accordingly, the momentum space 
integral 



FIG. 1: Enhancement of the particle density in the dynam- 
ically assisted Schwinger mechanism. The parameters of the 
adiabatic pulse are given by E\ = 0.1E C and wi ~ m/100. 
The different curves correspond to different values of E2 and 
we change 7 C or, equivalently, u>2- 



where the dot denotes the derivative with respect to t. 
Here, G(q,t) and H(q,t) are auxiliary functions and the 
remaining quantities are defined according to: 



UJ 2 {q,t) 



+ [q-eA(t)} 2 , W(q,t) 



eE(t)e ± 



(4) 



The initial conditions in the asymptotic past are given 
by F(q, -T) = G(q, -T) = H(q, -T) = 0. Note that we 
neglect the back reaction effect of created pairs on the 
electric field which is expected to be a good approxima- 
tion in the subcritical field strength regime. 



n[F) 



F(q,T)dq 



(2) 



can be interpreted as the density of created electron- 
positron pairs for qj_ = 0. 

The time evolution of the system is governed by an 
integro-differential equation [261 ] or, from a numerical 
perspective more conveniently, an ordinary differential 
equation system [27j |: 

F(q,t)=W(q,t)G(q,t) , (3a) 
G(q,t)=W(q,t)[l-F(q,t)]-2u;(q,t)H(q,t) , (3b) 
H(q,t) =2u{q,t)G{q, t) , (3c) 



1 Here a remark with respect to the recent investigation presented 
in Ref. [2G] is in order. We want to point out that presenting mo- 
mentum spectra at peak field values might be misleading as the 
distribution function can be related to the particle yield only at 
asymptotic times and vanishing fields (as has been acknowledged 
by the authors themselves) . In addition, based on our own simu- 
lations, we suspect that some of their results suffer from problems 
due to a lack of accuracy in the numerical computation. 



B. Dynamically assisted Schwinger mechanism 

Based on the quantum kinetic formalism, wc arc able 
to calculate the momentum spectrum as well as the num- 
ber of created electron-positron pairs for arbitrary field 
configurations. As mentioned previously, one type of field 
configuration which has attracted much interest recently 
consists of a superposition of a strong adiabatic pulse 
(index 1) with a weak anti- adiabatic pulse (index 2): 



E(t) = Ex sech 2 (wi£) + E 2 sech 2 (w 2 £) 



(5) 



The adiabaticity of a single pulse is quantified in terms 
of the Kcldysh parameter: 
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mu> 
~eE 



(0) 



This parameter discriminates between two different 
regimes: The pair creation process shows the charac- 
teristics of the Schwinger effect for 7 < 1 whereas the 
multiphoton regime is associated with 7 > 1. For the 
field configuration ([5]) we assume 71 <C 1 and 72 ^>> 1 or, 
equivalently, E c 3> E\ 3> E2 and uj\ <C W2 *C rn. It has 
been pointed out, however, that neither ji nor 72 is the 
relevant parameter when dealing with the dynamically 
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FIG. 2: Particle density for of a comb of 10 single pulses with 
parameter \Ei\ = 0.02_B C , uji = m/6 for the equal-sign configu- 
ration (solid) and the alternating-sign configuration (dashed). 
The particle number changes quasi-periodically by orders of 
magnitude as function of the inter-pulse time lag to- 



assisted Schwinger effect but rather a combined Keldysh 
parameter, which is composed of the dominant scales Q: 



7c 



mw2 



(7) 



In fact, the non-trivial interplay between different scales 
results in a subtle change of the pair production charac- 
teristics for properly chosen parameters. Most notably, 
it has been observed that F\ +2 (q,T) shows strong devia- 
tions compared to T) + F 2 (q, T) Consequently, 
this results in a strong enhancement of the particle yield, 
i. e. 71^1+2] can become orders of magnitude larger than 
n[F{\ +n[F 2 }. 

To illustrate this, we superimpose an adiabatic pulse 
with parameters E\ = 0.1 E c and uj\ <~ to/ 100 with an 
anti-adiabatic pulse for different values of E 2 as function 
of 7 C in FIG. [TJ It can be seen that the enhancement 
peaks around 7 C = as well as -f c ~ 2. The case j c = 
corresponds to a trivial effect as the second pulse turns 
from the anti-adiabatic to the adiabatic regime. Conse- 
quently, particle creation becomes enhanced simply due 
to a higher peak field strength Ei+E 2 . The peak around 
7 C ~ 2, however, is a nontrivial result of the analysis of 
the dynamically assisted Schwinger mechanism. Finally, 
the enhancement approaches unity for large values of j c . 
In this case, the second pulse completely dominates as it 
approaches the perturbativc pair creation threshold. 

This simple example nicely illustrates the general idea 
of the dynamically assisted Schwinger mechanism. How- 
ever, it is already rather challenging to investigate the 
whole four-dimensional parameter space {E\, E 2i uji,uj 2 }. 
This becomes even worse when we consider a more com- 
plicated field configuration such as a superposition of an 
arbitrary number of single pulses: 



E(t) = J2 E i sech a (w«[* + *o,il) 



(8) 



In this case, all pulses may have a different field strength 



Ei, frequency coi and time lag tg,,. Accordingly, a sys- 
tematic investigation becomes absolutely impracticable. 

In order to get an impression of the complicated struc- 
ture of the parameter space, we employ the idea pre- 
sented in Ref. [l8[ where it has been shown that a se- 
quence of alternating-sign pulses produces a Ramsey in- 
terferometer. Based on this observation, we investigate 
the alternating-sign configuration as well as the equal- 
sign configuration in FIG. [2] To be specific, we show the 
number of created particles for a comb of 10 single pulses 
with parameters \Ei\ = 0.02E C , = m/6 and to,i = i-to 
as function of the inter- pulse time lag to- Note that 
we have sign [£4] = +1 for the equal-sign configuration 
whereas we have sign[£*j] = ±1 in the alternating-sign 
configuration for i being even and odd, respectively. 

We observe some remarkable features upon varying the 
inter-pulsc time lag tQ. On the one hand, the number of 
created particles changes quasi-periodically by orders of 
magnitude for both the alternating-sign and the equal- 
sign configuration. This meets the expectations for the 
alternating-sign configuration according to the discussion 
in Rcf. [18j . however, the appearance of a pronounced in- 
terference effect for the equal-sign configuration has not 
been anticipated. On the other hand, the local extrema 
in the particle density show an out-of-phase behavior for 
these two configurations. We must admit at this point 
that we have not fully understood all effects for the equal- 
sign configuration yet. Nonetheless, the calculation in- 
dicates again that rather minor parameter changes can 
result in a tremendous enhancement of the particle yield. 
A more detailed investigation and discussion of the dif- 
ferences between the alternating-sign and equal-sign con- 
figurations can be found in Ref. [19( . 



III. PULSE SHAPING 

We have seen in the previous section that a proper 
parameter choice is crucial for maximizing the particle 
yield. However, we have also observed that a systematic 
scan of the parameter space becomes rather impractica- 
ble for field configurations more complicated than ([5|). In 
order to overcome this drawback, we will now introduce 
optimal control theory to increase the particle yield in a 
more systematic way. 



A. Optimal Control 

The task of maximizing the number of produced 
electron-positron pairs by varying the field configuration 
can be formulated as an optimization problem with con- 
straints. Moreover, it is obvious from a physical point 
of view that it is not sufficient to maximize n[F] without 
any restrictions. Due to the fact that the Schwinger effect 
is exponentially suppressed, for instance, an increase of 
the particle yield could be trivially achieved by increas- 
ing the field strength. As a consequence, we have to put 
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additional limitations on the admissible field configura- 
tions. Accordingly, the objective functional, which is to 
be minimized, is taken to be: 

J[F,A] = - 7 n[F]+f[A] . (9) 

Here 7 is some constant to be specified, A is the vec- 
tor potential and the functional f[A] has been added to 
restrict the set of admissible vector potentials A. The ob- 
jective functional is subject to the constraints ep,G,H(q, t) 
corresponding to the equations of motion ([3]): 

= e F = F-WG , (10a) 
= e G = G + WF + 2luH - W , (10b) 
= e H = H- 2loG . (10c) 

Apart from the inhomogeneous contribution W(q,t) to 
ec(q,t), these equations are determined by a skew- 
symmetric matrix which conserves the norm of the vector 
(F, G, H) in time. These constraints can be accounted for 
by introducing Lagrange multipliers fj,F,G,n(q,t) leading 
to the Lagrangian 

L(F,G,H,A,ii F ,lia,fiH)= (H) 
J[F, A] + (fi F , e F )a + (fia, e G )o + {hh, &h)q, , 



where (•, -)o denotes the L 2 inner product on 51 = 
R x [—T,T]. A vanishing gradient of the Lagrangian is 
a necessary condition for F, G, H, A to be a local min- 
imizer of J[F, A]. Setting the gradient with respect 
to (fip, A*G?> equal to zero reproduces the constraint 
equations ([ID]). 

More interestingly, the gradient with respect to 
(F, G, H) leads to the adjoint equations: 



= fip - W» G , (12a) 
= fi G + W[i F + 2u(j, H , (12b) 
= fi H - 2ufi G , (12c) 



with final conditions Hp(q,T) = 7 and /j,g(q,F) = 
HH{q, T) = 0. Apart from the inhomogeneous term these 
equations of motion for the Lagrange multipliers are iden- 
tical to ([TO]). 

Finally, the variation of the Lagrangian with respect to 
A can be calculated straightforwardly as well and gives, 
apart from time-independent boundary terms: 



SL 
5A 



5A 



2 ( PH G - n G H) 



q — eA e± d 1 

h II WgF - HfG - no) > dq 

oj ui dt 



(13) 



Assuming that the equations of motion ([10]) can be solved 
uniquely for any admissible potential A, it is useful to 
define the reduced cost functional 



J[A] = J[F(A),A] 



(14) 



With the aid of J [A] the formerly constrained optimiza- 
tion problem can be reformulated as an optimization 
problem without constraints. As a matter of fact, most 
optimization methods require at least knowledge about 
the gradient of the corresponding objective functional. In 
this respect, it can be shown that the gradient of the re- 
duced cost functional with respect to the control function 
A is given by 



6J[A] 
SA 



SL 
5A 



(15) 



if F,G,H is the unique solution of (|10| and fip , \ig^h 
fulfill the adjoint equations (fT2"]) . Accordingly, apart from 
solving the equations of motion forwards in time it is 
also necessary to evolve the adjoint equations backwards 
in time in order to calculate search directions for the 
reduced cost functional. 



B. Constraints 

The previous derivation has been general in the sense 
that we have not specified the functional f[A] in order 
to implement physical restrictions on the admissible field 
configurations. As a matter of fact, the possibilities are 
virtually unlimited but for the sake of simplicity we will 
discuss only two restrictions at this point. 

The first limitation concerns the electric field strength 
as discussed previously. To ensure boundedness of the 
electric field E — — A at any time 

E min (t) < -A(t) < E max (t) , (16) 

we implement the cost functional: 

fi[A] = -Jih[-A, -E min ] - l2 h[A, E max ] , (17) 

with 

T 

h[A,E} = J log (E(t) + A(t)j dt . (18) 

-T 
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Here 71,2 are positive numbers and both E m { n and E max 
should vanish at asymptotic times ±T to ensure an 
asymptotic particle interpretation in this limit. The func- 
tional will diverge fi[A] — > +00 if cither — A \ E min 
or —A /• E max on some subinterval of [—T,T]. Ac- 
cordingly, field configurations out of the bounds (fTr?|) 
are prevented as they would render the minimization of 
the objective functional impossible. Apart from time- 
independent boundary terms, its variation with respect 
to the control function is given by: 



5A 



-71 



A + E n 



(A + E n 
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A + E n 



(A + E n 



(19) 



Another sensible restriction might be put on the total 
amount of energy <f max . This constraint is again imple- 
mented via a functional 



with 



MA] = -73l0g(# m ax-^]) 



m = 2 J A\t)dt . 

-T 



(20) 



(21) 



Again, this functional will diverge f2[A] — > +00 if <? ma 
is approached. Consequently, the derivative is given by 

A 



SA 



-73 



<S<< 



S[A\ 



(22) 



C. Parametrized Field Configuration 

In certain situations it might be convenient to formu- 
late the optimization problem also on a finite dimensional 
space. As an example consider the field configuration © 
where we fix the functional form but still allow certain 
parameters to change. 

In order to account for such problems we assume that 
the vector potential is given by a function 



A{t) -> A(t;p) , pe 



(23) 



where p is the set of parameters which are to be opti- 
mized. For the field configuration ((5]), for instance, p 
would be the set of field strengths Ei , frequencies and 
time lags toj. The Lagrangian is then again defined ac- 
cording to 



L(F,G, H,p, fi F , fj, G , n H ) = 

J[F, A(p)] + (n F ,e F )n + (mg, e G )o 



(24) 



such that the equations of motion (|10p and adjoint equa- 
tions (Ti"2"|) are not changed apart from the replacement 
A(t) -> A(t;p). 

Moreover, the reduced cost functional becomes a func- 
tion of p 

J(p) =^[F(A(p)),A(p)] . (25) 
Its gradient V p J, which becomes a vector in R™, can 
again be obtained from the Lagrangian by evaluating its 
derivative with respect to p 
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V P J = V P L = V p /L4(p)] + e J J i 2 (p H G - n G H) ^— ^ V p A - e± (p G F - f i F G - MG ) V p (-^) I dtdq . (26) 



-00 -T 



for F,G,H and fj, F , fi G , fin being the solutions of the 
equations of motion (fT0|) and the adjoint equations (TTSj) . 
respectively. 

Further restrictions on the admissible values of p are 
again specified in terms of /[A(p)]. If we were to use the 
constraints ([16]) or (|20|) . respectively, the corresponding 
gradients would read 

and 

T 

v ' h =^hmI Av - Adt - ,28) 

-T 



D. Optimization 

Given the gradient of the reduced cost functional we 
are in the position to employ standard optimization 
methods. The general strategy is to approach a local 
minimizer of the reduced cost functional iteratively via 

Pk+i = Pk + a k d k , (29) 

where p k are the expansion coefficients of the vector po- 
tential, dk is a search direction and a k the step size in the 
k-th iteration step. In the infinite dimensional formula- 
tion, pk corresponds to the values of the vector potential 
at intermediate times A(t) whereas it represents the vec- 
tor of parameters p in the finite dimensional formulation. 

The search direction d k is taken such that the reduced 
cost functional is improved upon iteration, i.e. driven 
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towards a local minimum. In general, this condition im- 
plies that its inner product with the gradient should be 
negative. As a matter of fact, there are various possible 
choices for the search direction dk such as taking the neg- 
ative gradient - termed the steepest descent method - or 
multiplying the gradient with the negative inverse Hes- 
sian - called Newton's method. Other possibilities which 
are usually faster than steepest descent but also more 
involved include conjugate gradient methods or quasi- 
Newton methods. 

Finally, the step size Uk is chosen such that the func- 
tion 4>y. (a) = J(pk + adk) fulfills the strong Wolfe condi- 
tions [28[ 

M<*k) < 0*(O)+cia fc ^ fc (O) , (30a) 
|&(o*) I < ca 1^(0)1 , (30b) 

with < c\ < C2 < 1 at a^. Here, the derivative of 
4>k{ot) is given by the inner product (VJ (p& + adk) , dk)- 
In the infinite dimensional formulation, an L 2 inner prod- 
uct on [— T, T] is taken and V represents the functional 
derivative whereas the usual Euclidean inner product on 
l n with V being the gradient is employed for the finite 
dimensional formulation. 

Having specified the search direction dk - in our sim- 
ulations we will employ the steepest descent as well as 
the Fletcher-Reeves conjugate gradient method - there 
remains the task of finding an appropriate step size otk 
which fulfills the strong Wolfe conditions (|30|) . Again, we 
will employ a simple algorithm here [28]: Our strategy 
is to increase the step size a until the interval [0,a] con- 
tains admissible steps. As soon as this is achieved, we 
use bisection on this interval. 



IV. RESULTS 

In the following we present some results which are 
based on the optimal control formalism. For that pur- 
pose we initialize the system at a simple field configura- 
tion which is then modified according to the optimization 
procedure. The equations of motion (| 10[) and the adjoint 
equations (|12|) are solved with the aid of a Dormand- 
Prince Runge-Kutta integrator of order 8(5,3). For sim- 
plicity we restrict ourselves to parametrized field config- 
uration of section IIII CI 

A. Necessity of constraints 

As simple example we start with the Sauter field con- 
figuration 

A in {t) = -^tanh(wt) , (31a) 
E- m (t) = E scch 2 (wt) . (31b) 

For the moment we focus on the optimization in the one- 
dimensional space p — Eq whereas we keep u> fixed. As 
initial values we specify Eq = 0.1i? cr and uj = ra/15. 



B(t) 




t[l/m] 



FIG. 3: Electric field configuration: The initial configuration 
(solid) is characterized by Eq = O.liJcr- After a few iteration 
steps the field strength Eq converges towards the field strength 
constraint -Bo, max = 0.15i5cr. The final configuration (dashed) 
is thus characterized by Eo — -Eo.max. 

Upon performing the optimization we observe that the 
gradient V p J is indeed directed towards larger values of 
Eo- However, we cannot find permissible step lengths 
afc consistent with the strong Wolfe conditions (foTJf as 
the objective functional J(p) is not bounded from below. 
Accordingly, trying to minimize the objective functional 
becomes equivalent to taking the limit Eq — > oo. 

In order to set up a convergent optimization proce- 
dure we need to include some constraint which prevents 
the field strength from growing without restraint. This 
can be achieved by specifying a functional /L4(p)] which 
has been introduced in @ in order to restrict the ad- 
missible vector potentials. Therefore, we start with the 
same initial configuration (|31| and employ, for instance, 
the field strength constraint (|16| with: 

Eo,max = — -Eo,min = 0.15i? cr . (32) 

As a consequence, the objective functional J(p) be- 
comes bounded from below. The resulting field configura- 
tion after a few iteration steps is displayed in FIG. [3l As 
intended, the maximum field strength becomes bounded 
and the optimization procedure converges. This shows 
that it is sometimes absolutely necessary to choose an 
appropriate constraint in order to obtain sensible opti- 
mization results. 

Note that we could also employ the energy constraint 
(f20"|) in order to prevent the electric field strength from in- 
creasing without bounds. This type of constraint is qual- 
itatively similar to a field strength constraint as choosing 
a maximum energy 4iax effectively determines a maxi- 
mal field strength as well. 

B. Optimization towards the threshold 

We consider again the Sauter field configuration ([3T]) 
but this time we fix Eq = 0.1i? cr . Accordingly, we 
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FIG. 4: Particle number n[F] for the Sauter field with 
Eo — 0.1-Ecr as a function of cj. The maximum near the 
perturbative threshold cj op t ~ 1. 83m is also predicted by our 
optimal control procedure. The arrows indicate the downhill 
directions of the objective functional J = — yn[F]. 

perform an optimization in the one-dimensional space 
p = uj. Note that we do not need an additional constraint 
in this case as the field strength is specified anyway. 

In order to perform the optimization, we choose rather 
different initial values uj = m/20 and uj = 10m, respec- 
tively. Upon optimization, we find that the optimal value 
converges towards w pt ~ 1.83m for both initial values. 
As a matter of fact, we may confirm that this value corre- 
sponds to the maximum in the particle number near the 
perturbative threshold by comparison with an analytic 
result. 

In fact, the asymptotic distribution F (q, T) for the 
Sauter field is given by [29| : 

= 2sinh(^[2 £ + / i-^])sinh(^[2 £ - / i + ^]) 
W ' ' sinh(^)sinh(^) 

(33) 

with e = eE /uj and 

M = V /e l + (9 + e ) 2 . " = V c i + (9- £ ) 2 ■ ( 34 ) 

Accordingly, we can explicitly evaluate n[F] which is dis- 
played in FIG. [U Most notably, this curve shows a max- 
imum near the perturbative threshold at uj op t ~ 1.83m 
which coincides with our numerical result. For larger 
values uj > w opt , the particle number decreases with a 
power law us~ 2 . On the other hand, for smaller values 
m/50 < uj < uj opt we observe the characteristic transition 
between the multiphoton regime to the Schwinger regime 
which is accompanied by an order-of-magnitude change 
in the particle numbers. Moreover, we find an additional 
maximum at the boundary uj — > which is approached 
with a power law lu^ 1 . In this limit the Sauter field (|3ip 
adiabatically changes into a static field E(t) = Eq such 
that n[F] diverges. 

Illustrated by this simple example we want to discuss 
a numerical challenge which may occur in more compli- 
cated situations as well: Starting with an initial value 



uj < m/50 one observes that the local downhill direc- 
tion of the objective functional J = —^n[F] is directed 
towards uj — > 0. This means that any attempt to start 
the optimization procedure in this regime should finally 
drive u — > instead of u — > uj opt . Note, however, that 
the objective functional J is rather flat in this regime 
so that the gradient V p J is close to numerical precision. 
Accordingly, the whole optimization procedure becomes 
rather untrustworthy. 

C. Local maxima vs. global maximum 

In order to discuss another general issue of optimal 
control theory, we turn to a more complicated field con- 
figuration than the Sauter field ([3Tj) . To this end we 
reconsider the field configuration (|8|) corresponding to a 
superposition of arbitrary single pulses. The need to in- 
vestigate this kind of field configuration more efficiently 
was in fact one of our main motivations for the applica- 
tion of the optimal control formalism. 

We now re-investigate the comb of 10 single pulses in 
both the equal-sign and the alternating-sign configura- 
tion: 

10 

Ain(t) = -^fftanh(o;[t-i-to]) , (35a) 

10 

E iD (t) = ]T E, sech 2 (w[i - i ■ t }) . (35b) 

i=l 

Note again sign[i?j] = +1 for the equal-sign configura- 
tion and sign[_Ey = ±1 in the alternating-sign configu- 
ration for i being even and odd, respectively. Based on 
our experience we expect the optimization to drive Ei 
towards higher field strengths and w towards the per- 
turbative threshold, respectively. Accordingly, we fix 
\Ei\ = 0.02i? cr and ui = m/6 such that the the optimiza- 
tion problem is again formulated in a one-dimensional 
space p = to °f the inter-pulse time lag. 

The equal-sign configuration is investigated in TAB. U 
whereas results for the alternating-sign configuration are 
shown in TAB. [TTJ We display the initial values £o,init, 
the corresponding optimized values to, opt as well as the 
optimal values fo.comp which arc deduced from FIG. [2] 
Most notably, we observe that the optimization proce- 
dure does not yield a unique optimal value corresponding 
to a global maximum of n[F] but a whole slew of different 



Initial to,init 


Optimized to, opt 


Comparison to,com P 


15 


15.86 


15.87 


17 


15.75 


15.87 


19 


18.96 


18.98 


21 


22.12 


22.10 



TABLE I: Optimization of the inter-pulse time lag to for the 
equal-sign configuration. The values to, op t are in good agree- 
ment with to,com P which are deduced from FIG. [5] 
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Initial to,init 


Optimized to, opt 


Comparison to.comp 


Initial (i?2,W2) 


Optimized (#2,^2) 


Optimized n[F] 


16 


15.43 


15.46 


(10 4 ,10) 


(0.0205, 1.830) 


1.176 • 10~ 4 


1 8 


1 7 no 

1 1 .UU 


1 7 1 Q 




(C\ 090/1 1 R r )7\ 


1 17c: in -4 

1.1(0 • 1U 


19 


18.50 


18.57 


(0.005, 10) 


(0.0205, 1.819) 


1.182 ■ 10~ 4 


20 


20.50 


20.30 


(0.005, 1/2) 


(0.0204, 1.827) 


1.175 • 10~ 4 


21 


21.75 


21.68 


(0.02, 10) 


(0.0205, 1.833) 


1.181 • 10~ 4 


23 


23.50 


23.43 


(0.02, 1/2) 


(0.0204, 1.823) 


1.175 • 10~ 4 



TABLE II: Optimization of the inter-pulse time lag to for the 
alternating-sign configuration. Note that we resolve both the 
dominant and the sub-dominant local maxima. 



optimal values corresponding to different local maxima of 
n[F]. The reason for this is the fact that optimal con- 
trol theory deals with the determination of local extrema 
whereas the identification of the global extremum is be- 
yond its scope. In order to search for a global maximum 
of n[F ] it is necessary to perform the optimization pro- 
cedure for various initial values. In this respect, a proper 
choice of the initial values can be most crucial: A poor 
choice may put us in the basin of attraction of a 'wrong' 
minimum of the objective functional J(p). 

Moreover, if we compare the optimal values io.opt 
with io.comp, we observe quite reasonable agreement even 
though we do not obtain full convergence. It can only 
be speculated why this happens, however, it seems that 
our algorithm tends to considerably slow down once we 
are close to the local maximum of n[F]: Starting from a 
given initial value we proceed in very few iteration steps 
towards this maximum, however, once we are close to 
it the gradient V p J becomes rather small. As a con- 
sequence, its numerical evaluation becomes error prone 
such that the optimization procedure is terminated even 
though full convergence has not yet been achieved. 



D. Multidimensional optimization 

We now turn to higher-dimensional parameter spaces 
after discussing only one-dimensional optimization prob- 
lems so far. In the one-dimensional case, the parameter 
pel could only be increased or decreased depending 
on the sign of the gradient V P J. Upon increasing the 
dimensionality of the parameter space R™, however, the 
gradient becomes an n-dimensional vector and, associ- 
ated with that, several subtleties arise. Most notably, we 
will discuss the issue of scaling which is crucial in order 
to obtain sensible optimization results. 

For simplicity, we restrict ourselves to a two- 
dimensional parameter space in the following. To be 
specific, we investigate the double pulse configuration: 



An(t) = -^tanh(wit) - ^tanh(w 2 t) 
E in (t) = Ex scch 2 (o;ii) + E 2 sech 2 (w 2 i) 



(36a) 
(36b) 



TABLE III: Two-dimensional optimization for p = (£2,^2) 
in units of E CI and m, respectively. We display different initial 
configurations, the resulting optimized configuration as well 
as the particle number n[F]. 



upon variation of E 2 such that we are concerned with a 
two-dimensional optimization problem for p = (i? 2 ,w 2 ). 
Note that in order to guarantee \E(t)\ < E max for all t, 
we need to include a field strength constraint p^|) as well. 

In order to perform the optimization, we choose rather 
different initial configurations as displayed in TAB. IIIII 
In this table one can also see that all these configura- 
tions converge quite well towards an optimal configura- 
tion p pt ~ (0.0203-E cr , 1.82m). One of the initial con- 
figurations as well as the optimized field configuration is 
displayed in Fig. [5l 

The behavior of the optimized field configuration is 
in accordance with our previous investigations of one- 
dimensional problems. We have seen that a maximum 
particle yield is obtained for frequencies which are close 
to the perturbative threshold. In accordance with that, 
the frequency ui 2 is driven towards this threshold u; 2 . op t ~ 
1.82m. Moreover, the field strength E2 is maximized 
subject to the field strength constraint (fl"6]) as well as 
the condition E1+E2 = E max . This means, on the other 
hand, that the field strength E\ < ultimately takes a 
negative value. 

We have already mentioned, that the optimization 
in higher dimensions is associated with several sub- 
tleties, most notably the issue of scaling [H}: Given an 



E(t) 




t[l/m] 



We keep Wi = m/10 fixed and demand E x + E 2 = £ max 
with -E max = 0.01i?cr- In this configuration, E\ changes 



FIG. 5: Electric field configuration: For simplicity we show 
only one of the initial configurations (solid) which is character- 
ized by pi n = (10 -4 -E cr , m/2). The optimized configuration 
(dashed) converges towards p opt ~ (0.203-E cr , 1.82m). 
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Iteration Steps 



FIG. 6: Rate of convergence of the particle number n[F] for 
different values of the diagonal scaling exponent x defined 
in (|37[) . In each case we start with the initial configuration 
Pin = (0.02_E cr , 10m). The connecting lines are included to 
guide the eye. 

optimization problem in higher dimensions it might hap- 
pen that the objective functional J(p) changes rather 
differently as a function of the various parameters Pi. 
Consequently, the various components of the gradient 
V p J(p) may also differ substantially so that the descent 
direction is dominated by a specific subset of parameters. 
If this is the case, the problem is said to be poorly scaled 
as the rate of convergence towards p opt might decrease 
significantly. As a matter of fact, algorithms such as the 
steepest descent method are very sensitive to poor scal- 
ing behavior. To remedy this shortcoming it is indicated 
to perform diagonal scaling, i. e. redefine p — > p' such 
that the objective functional J(p') is better balanced as 
a function of p'. 

As a matter of fact, our objective functional (j^J turns 
out to be poorly scaled as a function of p = (i^Wa)- 
Accordingly, we need to perform diagonal scaling in order 
to optimize the convergence properties of our algorithm: 

lo' 2 = 10- x lu 2 . (37) 

In Fig.[6]we investigate the rate of convergence for a given 
initial configuration pi n = (0.02E CI , 10m) for different 
values of the diagonal scaling exponent x. We find that 
the rate of convergence, i. e. the number of iteration steps 
which have to be taken until n[F] converges towards its 
maximum value, strongly depends on the choice of x. 

As a matter of fact, for x = the ^-component of 
the gradient V p J(p) is orders of magnitude larger than 
its o^-componcnt. Accordingly, in the course of the op- 
timization the field strength quickly approaches its op- 
timal value £2, opt, however, the optimization of u>2 does 
not make much progress. Due to numerical inaccuracies 
the ^-component of the gradient V p J(p) remains larger 
than its W2-component even if we are close to i?2,opt- 
Hence, we are finally stuck at a point in parameter space 
at which n[F] is still way below its optimal value. 

On the other hand, we observe that the rate of conver- 
gence improves significantly upon increasing the diagonal 



scaling exponent x. This clearly shows that the issue of 
scaling is a most crucial one when we want to perform op- 
timization in higher dimensions. Most importantly, poor 
scaling might even prevent the algorithm from converging 
towards the optimal value. 



V. CONCLUSION 

Based on optimal control theory, we performed a sys- 
tematic pulse shaping analysis for electron-positron pair 
creation in strong external fields. To this end, we de- 
rived the optimal control equations for both finite and 
infinite dimensional parameter spaces in the framework 
of quantum kinetic theory. As this was the first study of 
this kind, we then focused on the finite dimensional for- 
mulation for simplicity. We convincingly demonstrated 
that optimal control theory provides a unique means for 
systematically maximizing the particle yield. 

The need for pulse shaping has been pointed out previ- 
ously and recent investigations have again indicated that 
rather minor parameter changes may increase the parti- 
cle yield by orders of magnitude [l8[ . Due to the fact that 
a systematic scan of the parameter space becomes rather 
impracticable soon, it is indicated to perform these pulse 
shaping analysis by means of an automated tool such as 
optimal control theory: Within this framework, the par- 
ticle yield is driven towards a maximum by an iterative 
variation of the field configuration. 

For the Sauter field configuration, where a comparison 
with an analytic result is possible [29| . we reproduced 
the trivial maximization directions: The particle yield 
is maximized by increasing the field strength at a given 
frequency or by driving the frequency towards the pertur- 
bative threshold at a given field strength. For more com- 
plicated configurations, however, several subtleties arose 
which will also have to be accounted for in the future: 

First, optimal control theory deals with the determina- 
tion of local maxima of the objective functional whereas 
the identification of the global maximum is beyond its 
scope. Accordingly, a poor choice of initial conditions 
may put us in the basin of attraction of a 'wrong' max- 
imum. As it is not possible to comb the whole potential 
landscape for all local maxima, we still rely on educated 
guesses or physical intuition for the initial conditions. 

Secondly, we have seen that a maximization of the par- 
ticle yield can always be achieved by increasing the field 
strength or choosing the frequency near the perturbative 
threshold. In order to exclude these trivial maximization 
directions, we need to apply appropriate constraints. We 
demonstrated the action of a field strength constraint 
convincingly whereas the implementation of a frequency 
constraint in the current formulation has not been at- 
tained yet. 

Finally, one of the most severe problems we encoun- 
tered during our numerical investigation was the poor 
scaling behavior of the objective functional. Due to 
the fact that we employed algorithms which arc rather 
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sensitive to poor scaling (steepest descent, conjugate gra- 
dient without preconditioning), it was essential to per- 
form diagonal scaling even in a two-dimensional param- 
eter space. As this problem becomes even more serious 
in higher dimensions it seems to be inevitable to replace 
these simple algorithms in the future by more advanced 
ones such as quasi-Newton method, which is not affected 
by the poor scaling behavior. 

In our opinion, the issue of pulse shaping might become 
a crucial one once the new generation of high-intensity 
laser facilities such as the Extreme Light Infrastructure 
(ELI) or the European XFEL are operating. Due to the 
fact that small parameter changes might enhance the par- 
ticle yield by orders of magnitude, it is surely worthwhile 
to scan the parameter space for optimal field configura- 
tions. In such optimized field configuration, the obser- 
vation of electron-positron pair creation might become 
possible even in the sub-critical field strengths regime. 

Due to the fact that the results presented in this pa- 
per have been based on rather simple parametrized field 
configurations, a natural next step would be to set up 



the optimization algorithm for more realistic field config- 
urations. As a matter of fact, first steps in this direction 
have already been taken [l9| . In the long run, these stud- 
ies should provide us with clear predictions on optimized 
field configurations for which the suppression of the pair 
production process becomes partially lifted. Given that 
such optimal field configurations exist, an experimen- 
tal verification of non-perturbative electron-positron pair 
creation seems realistic within the next decades. 



Acknowledgments 

MM and CK were funded by the Austrian Science 
Fund, FWF, through the Doctoral Program on Hadrons 
in Vacuum, Nuclei, and Stars (FWF DK W1203-N16). 
FH was supported by the Alexander-von-Humboldt 
Foundation. We thank the research core area "Model- 
ing and Simulation" for support. 



F. Sauter, Z. Phys. 69, 742 (1931). [18 
W. Heisenberg and H. Euler, Z. Phys. 98, 714 (1936). 
J. S. Schwinger, Phys. Rev. 82, 664 (1951). [19 
A. Ringwald, Phys. Lett. B510, 107 (2001). 

E. Brezin and C. Itzykson, Phys. Rev. D2, 1191 (1970). [20 
V. Popov, Sov. Phys. JETP 34, 709 (1972). 

R. Alkofer, M. Hecht, C. D. Roberts, S. Schmidt, and [21 

D. Vinnik, Phys. Rev. Lett. 87, 193902 (2001). 

C. Bamber, S. Boege, T. Koffas, T. Kotseroglou, [22 

A. Melissinos, et al., Phys. Rev. D60, 092004 (1999). 

R. Schutzhold, H. Gies, and G. Dunne, Phys. Rev. Lett. [23 

101, 130404 (2008). 

A. Di Piazza, E. Lotstedt, A. Milstein, and C. Keitel, [24 
Phys. Rev. Lett. 103, 170403 (2009). 

G. V. Dunne, H. Gies, and R. Schutzhold, Phys. Rev. [25 
D80, 111301 (2009). 

A. Monin and M. Voloshin, Phys. Rev. D81, 025001 
(2010). [26 
S. Bulanov, V. Mur, N. Narozhny, J. Nees, and V. Popov, 
Phys. Rev. Lett. 104, 220404 (2010). [27 
M. Orthaber, F. Hebenstreit, and R. Alkofer, Phys. Lett. 
B698, 80 (2011). [28 
C. Fey and R. Schutzhold, Phys. Rev. D85, 025004 
(2012). [29 

F. Hebenstreit, R. Alkofer, G. V. Dunne, and H. Gies, 
Phys. Rev. Lett. 102, 150404 (2009). 
C. K. Dumlu, Phys. Rev. D82, 045007 (2010). 



E. Akkermans and G. V. Dunne, Phys. Rev. Lett. 108, 
030401 (2012). 

C. Kohlfurst, arXiv: 1212.0880 [hep-ph] (Master's thesis, 
Graz University, 2012). 

A. Nuriman, B.-S. Xie, Z.-L. Li, and D. Sayipjamal, 
Phys. Lett. B717, 465 (2012). 

A. P. Peirce, M. A. Dahleh, and H. Rabitz, Phys. Rev. 
A37, 4950 (1988). 

J. P. Palao and R. Kosloff, Phys. Rev. Lett. 89, 188301 
(2002). 

J. Grond, G. von Winckel, J. Schmiedmayer, and U. Ho- 
henester, Phys. Rev. A80, 053625 (2009). 

F. Hebenstreit, R. Alkofer, and H. Gies, Phys. Rev. Lett. 
107, 180403 (2011). 

D. Blaschke, B. Kaempfer, A. Panferov, A. Prozorkevich, 
and S. Smolyansky, arXiv: 1205.3154 [physics. plasm-ph] 
(2012). 

S. Schmidt, D. Blaschke, G. Ropke, S. Smolyansky, 
A. Prozorkevich, et al., Int. J. Mod. Phys. E7, 709 (1998). 
J. C. Bloch, V. Mizerny, A. Prozorkevich, C. D. Roberts, 
S. Schmidt, et al., Phys. Rev. D60, 116011 (1999). 
J. Nocedal and S. J. Wright, Numerical Optimization 
(Springer, 1999). 

F. Hebenstreit, R. Alkofer, and H. Gies, Phys. Rev. D82, 
105026 (2010). 



